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文章 编号 :1005-3085(2010)04-0663-06 
非 线 性 常 微分 方程 边 值 问题 的 最 优化 算法 * 


REI, RO M, 吴 海 涛 
(沈阳 建筑 大 学 理学 院 ， 沈 阳 110168) 
摘 要 : 本 文 提 出 了 一 种 新 的 优化 方法 以 精确 地 解决 微分 方程 的 边 值 问题 。 针 对 非 线 性 微分 方程 边 值 问题 
打靶 法 的 不 足 ， 以 未 知 的 部 分 初始 条 件 为 设计 变量 ， 以 给 定 的 边界 函数 值 与 设计 变量 之 间 的 关 
系 ， 形 成 复杂 和 岗 套 式 程式 化 的 目标 函数 ， 建 立 了 求解 未 知 初始 条 件 的 优化 问题 ， 并 采用 Visual 
Basic 6.0 语 言 编制 了 求解 未 知 初始 条 件 的 程序 ， 给 出 了 典型 算 例 ， 并 通过 比较 ， 验 证 了 本 文 求解 


方法 的 正确 性 。 
关键 词 : 非 线 性 微分 方程 ， 边 值 问 题 ， 优 化 方法 ， 程 序 设 计 
分 类 号 : AMS(2000) 34B15 中 图 分 类 号 : 0175.8; O224 文献 标识 码 : A 
1 引言 


非 线性 常 微分 方程 边 值 问题 是 工程 中 常见 的 理论 和 实际 问题 [" 引 ， 通 常 采用 有 限 差分 法 和 打 
惑 法 分 析 求 解 。 目 前 应 用 最 多 的 求解 方法 是 打靶 法 申 。 无 论 单 级 和 多 级 打靶 求解 方法 ， 通 过 下 
面 步骤 求解 : 1) 选取 待 求 的 初 值 ，2) 采用 微分 方程 初 值 问 题 算法 按 步 迭代 计算 ， 获 得 终点 边 
界 值 ，3) 应 用 计算 边界 值 与 给 定 真 实 边界 值 的 误差 函数 ， 重 新 调整 待 求 的 初 值 ， 直 到 使 得 误差 
函数 达到 一 定 精度 为 止 。 总 的 看 这 种 计算 过 程 比较 繁琐 。 最 优化 方法 是 能 良好 解决 工程 实际 问 
题 的 数学 方法 中。 笔者 曾 提出 动态 设计 变量 优化 方法 通过 巧妙 进行 动态 设计 变量 处 理 ， 构 造 目 
标 函 数 ， 已 经 解决 了 许多 工程 实际 问题 Bi0 。 本 文 针 对 非 线性 微分 方程 边 值 问题 ， 提 出 一 种 以 
待 求 初 值 为 设计 变量 ， 以 误差 函数 构造 目标 函数 的 直接 最 优化 算法 。 并 通过 编制 最 优化 算法 通 
用 程序 实现 非 线 性 常 微分 方程 边 值 问题 的 精确 求解 。 


2 ”微分 方程 边 值 问题 的 最 优化 算法 原理 与 程序 设计 


2.1 微分 方程 边 值 问 题 
设 任意 一 个 nn 阶 非 线 性 常 微分 方程 


y™ -— f(zywy.w D^ 9070. (1) 
边界 条 件 
y? (a) = y(9, k=ii,i2,.… Ini y (b) = y, l= ji j2: jm Ry (2) 


其 中 有 ,i2,… sini AI ji, jz iima Jna 为 小 于 nn 的 互 异 正 整 数 ， Hn Tn25-—n. 
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求解 方程 (1) 边界 问题 的 实质 是 寻找 一 种 可 行 的 方法 ， 确 定 其 中 一 个 边界 点 的 所 有 y, y, y, 
…- Vy 07D 的 值 , 使 得 问题 成 为 初 值 问题 .对 于 z = a 点 , BAM y (a) = y9, k= inin in 
ni ^B, ifüyO (a), l= 41,32, dms ,jns，mn2 个 值 是 未 知 的 。 边 值 问 题 算法 将 确定 yD (a), 
l2 ju jaime ,jns， 这 nz 个 未 知 量 。 下 面 将 给 出 求法 与 算 例 。 

2.0 算法 原理 

Hs Wu-yy-—ycyam Y"-)， 将 (1) 转 化 为 状态 方程 


ýı ys, 
Y2 = Y3, 
(3) 
Yn-1 = Vn; 
Ün > f(z, Y1: 2， n1); 
边界 条 件 转化 为 
yx(a) = Yk, k = l,i2,.， sinis y (b) = yi, l= ji j2 ey 
构建 最 优化 问题 
min f(z), (4) 
f(z) = 》 [us (5,9) — Yim (0)] (5) 
m-—1 


设计 变量 z = [z1, z2,… ,zns] 1， 表 示 待 求 的 初始 值 y(a), 1 = ju jac jazo yj (0 Mm = 
1,2,.… ,mn2， 为 含有 待 求 设计 变量 的 复杂 妖 套 显 式 的 函数 。 其 计算 借助 于 数值 积分 Runge- 
Kuttal! 算法 ， 由 程式 按 步 迭代 计算 获得 的 边界 z = b 处 相应 各 阶 函 数值 。 这 样 ， 该 
问题 转化 为 寻求 待 求 设计 变量 ， 使 目标 函数 获得 极 小 值 问题 。 其 理想 解 条件 为 f(z) 一 
0, Hüy(z,b) = w(b), L= 二,j2，,… ,jnsa。 因 此 其 优化 求解 过 程 相当 于 求解 非 线 性 方程 
组 yj (z, b) = yj, (b), M= 1,2,… ,n2 的 解 。 

2.3 ”目标 函数 的 形成 程序 

目标 函数 程式 化 的 形成 步 又 为 :给 定 自 变 量 z， 计 算 步 长 为 h。 给 Tz € [ab] 4 EE. h= 
(b —a)/M, zy =a +kh, k=0,1,2,...,M, zo =a, zy = b. 

1) 由 已 知 初始 条 件 yx(a) = vx(xo) = ye k = 1,2,… in, 和 待 求 初始 条 件 yj,, (x0) = 
za, m 二 1,2,.… in: 采用 R-K 方 法 进行 一 步 计算 ， 获 得 y(71), k— 1,2,… ,n， 它 们 是 关 
于 z 的 函数 ; 

2) Hiyk(zi), k=1,2, ,n. WH yk(z2), k=1,2, ,n， 它 们 是 关于 zz 的 函数 ; 

3) Hiyk(z) k-12,--,n, WS y(zai) 上 二 1,2,… ,n， 它 们 是 关于 zz 的 函数 ; 

4) Biw(ru-i) k=1,2, n, 计算 yr (xm), k—-1,2,---,m 

5) 最 后 可 以 形成 目标 函数 


"n2 


f(z) = Y (ui, (2,9) — y. (9). 


i=1 


2.4 ”坐标 轮换 方法 概述 
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下 面 给 出 以 二 维 问题 为 例 的 多 维 坐 标 轮换 法 优化 过 程 计算 原理 ， 参 见 图 1， 在 x1Ozs 平 
面 上 ， 任 取 初 始点 z(@) = [O OT, MOHR, ARMIER, r = r 固定 
不 变 ， 仅 改变 zi 使 目标 函数 值 下 降 ， 由 min f(z1, zx) 一 维 搜索 获得 相应 的 最 优点 2 = 
[xD aT. RERA n (P 为 起 点 , z1 aU 固定 不 变 , 仅 改 变 rz, 沿 zz 轴 方 向 ， 由 min. f (a, 
23) 一 维 搜索 获得 相应 的 最 优点 z9) = [zt0D,zt017， 完 成 一 轮 搜索 ， 进 行 收 敛 条 件 判定 ， 
当 满 足 精度 时 ， 停 止 欠 代 ， 和 否则 从 中 出 发 进行 下 轮 搜索 计算 ， 直 到 满足 精度 获得 最 优 
点 zf = [z+,z3]7 为止 。 对 于 nn 维 问题 有 同样 的 计算 过 程 。 





图 1: 梯度 法 搜索 过 程 示意 图 


2.5 程序 组 成 

论文 采用 Visual Basic 编制 计算 程序 [1 ， 程 序 模块 组 成 含有 : 1) 主 程序 ，2) 一 维 搜索 子 程 
序 段 ，3) 坐标 轮换 法 子 段 ，4) 动态 设计 变量 目标 函数 子 段 。 

2.6 ” 解 的 存在 性 与 精度 问题 

1) 解 的 存在 性 问题 ， 本 文 的 优化 方法 与 差分 法 和 打靶 法 相 比 ， 仍 属 数值 算法 范畴 ， 是 基 
于 微分 方程 边界 问题 的 解 存在 为 假设 前 提 条 件 下 的 一 种 更 为 精确 算法 。 本 文 不 讨论 微分 方程 边 
界 问题 解 不 存在 的 情况 。 

2) 步 长 对 精度 的 影响 : 本 文 微 分 方程 迭代 计算 过 程 ， 均 采用 四 阶 显 式 Runge-Kutta 法 ， 
它 的 稳定 区 域 为 (-2.78,0)， 通 常 步 长 越 小 越 稳定 。 在 足够 小 步 长 条 件 下 ， 计 算 误差 随 步 骤 增 
加 而 不 断 减 小 ， 具 有 稳定 与 相 容 性 。 将 保证 优化 计算 中 获得 的 解 具有 收敛 性 。 


3 AHH 


算 例 1 线性 常 微分 方程 边 值 问 题 y +y=0, y(0) = 0, y(7/2)= 1. 

本 题 为 非 线 性 常 微分 方程 一 个 特例 。 设 yl = y, ya = y ， 则 方程 转化 为 = yo, yo = 
一 W。 该 题 设计 变量 维 数 n = 1， 待 求 变量 zi = y (0) = yo(0)， 取 步 长 h = r/2/10000， 初 
值 z1(0) = 0.5。 设 定 黄金 分 割 法 精度 为 e = 10-， 通 过 程序 优化 计算 zi = 0.99997。 此 时 目标 
函数 精度 达到 e = 10-10。 本 题 的 真实 解 y(z) = sin x。 准确 值 为 y (0) = cos0 = 1.00000， 优 化 
结果 表明 具有 较 高 的 计算 精度 。 表 1 列 出 选取 不 同 迭 代步 长 hh 时 ， 应 用 程序 优化 计算 的 设计 变 
量 结果 和 误差 情况 。 可 见 和 迭代 步 长 越 小 ， 将 会 导致 优化 计算 结果 的 准确 程度 越 高 。 
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表 1: 算 例 1 中 不 同步 长 对 于 优化 结果 的 影响 










误差 百分比 优化 计算 的 z ”误差 百分比 
7/2/100 0.99658 


71/2/10 0.96635 3.365 


优化 计算 的 对 
7/2/10000 0.99997 
1/2/1000 0.99965 0.035 












算 例 2 含 一 个 未 知 量 的 非 线性 常 微分 方程 边 值 问 题 
y —Qü-z27)y-1, y(0-1, y(1)=3. 


Wy y =y; UI EAE y, — yo, yo = 1- (1 22)yi- 

本 题 设计 变量 维 数 为 n = 1， 此 时 z1 = y (0) = yo(0)。 取 步 长 为 h = 0.0001， 设 目标 函数 
精度 为 e = 10-5， 初 值 z1(0) = 0.5， 程 序 优 化 计算 得 z1 = 0.6419， 即 y (0) = 0.6419。 以 此 初 
值 ， 应 用 R-K 方 法 计算 ， 列 出 z 间隔 为 0.01-0.05 的 精确 解 ， 见 表 2。 图 2 为 优化 计算 后 初 值 所 
确定 的 函数 和 导数 随 坐 标 变化 过 程 。 


表 2: 算 例 2 中 应 用 优化 后 初 值 所 计算 的 浮 数 值 











2.7284 
0.94 2.7653 3.7184 
0.95 2.8028 3.7811 
0.96 2.8409 3.8451 
0.97 2.8797 3.9103 
0.98 2.9191 3.9769 
0.99 2.9592 4.0448 
1.00 3.0000 4.1141 


0.10 1.0743 0.8458 





1 0.6419 
0.01 1.0065 0.6619 
0.02 1.0132 0.682 

0.03 1.0202 0.7022 
0.04 1.0273 0.7224 
0.05 1.0346 0.7428 
0.06 1.0421 0.7632 
0.07 1.0499 0.7837 
0.08 1.0578 0.8043 
1.066 0.825 









0.9515 
1.0604 
1.1733 
1.2909 
1.4139 
1.5432 
1.6799 
1.8250 





0.15 
0.20 
0.25 
0.30 
0.35 
0.40 
0.45 
0.50 





1.1192 
1.1695 
1.2253 
1.2869 
1.3545 
1.4284 
1.5089 
1.5965 
























































E e 


H2: 优化 计算 后 初 值 所 确定 的 函数 和 导数 随 坐 标的 变化 过 程 
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40/3 ”Duffing 方 程 y 十 y 十 0.5y3 = 0. 

首先 以 初始 条 件 y(0) = 1, y (0) = 2， 按 步 长 为 h = 0.0001， 以 R-K 算 法 直接 计算 出 x = 
1 时 ，y(1) = 1.4524, y (1) = 一 1.4243。 形 成 含 二 个 未 知 量 的 非 线 性 常 微分 方程 边 值 问 题 yy + 
y t 0.5y3 = 0， 边 值 条 件 y(1) = 1.4524, y (1) = 一 1.4243， 用 优化 算法 反 求 y(0) 和 wy (0)。 本 
题 设计 变量 为 n = 2。 设 y(0) = z, y (0) = z2， 步 长 仍 取 为 h = 0.0001, 4ER 29 = 0.5, 29 = 
0.5， 黄 金 分 割 法 精度 为 e = 10-5， 多 维 目标 函数 精度 为 e = 10-10。 应 用 优化 程序 分 析 ， 经 
过 17 轮 x2 = 34 次 优化 计算 ， 目 标 函 数 达到 设 定 精度 。 得 z = 0.9999, 2; = 1.9990， 结 果 对 
真实 初 值 具有 非常 高 的 逼近 精度 ， 表 明 方 法 有 效 。 将 优化 过 程 中 设计 变量 和 目标 函数 的 变化 过 
程 列 在 表 3。 


表 3: 和 工 例 3 优 化 求解 过 程 





优化 次 数 zı 22 opf 优化 次 数 zu z2 opf 












0 0.5 0.5 2.0134656620 18 1.0012 1.9963 0.0000030899 
2 1.2116 1.6028 . 0.0963231334 20 1.0006 | 1.9975 . 0.0000009047 
4 1.0044 1.8133 . 0.0174709276 22 1.0003 1.9982  0.0000003382 
6 1.0475 . 1.903 . 0.0042758058 24 1.0001 1.9986  0.0000000665 
8 1.0254 1.9474 0.0011996556 26 1.0000 1.9987 . 0.0000000259 


28 1.0000 | 1.9988 . 0.0000000090 
30 1.0000 1.9989 . 0.0000000032 
32 0.9999 1.9989 . 0.0000000011 
34 0.9999 . 1.9990  0.0000000001 


10 1.0138 1.9705 . 0.0003526104 
12 1.0077 1.9835  0.0001110789 
14 1.0042 1.9903 . 0.0000337805 
16 1.0023 1.9942 . 0.0000101650 





4 结论 


本 文 基于 优化 原理 ， 提 出 了 非 线性 常 微分 方程 边 值 问题 的 直接 优化 算法 。 编 制 边 值 问题 的 
通用 的 求解 程序 。 真 实地 分 析 求 解 了 待 求 未 知 量 为 1 个 和 2 个 的 线性 和 非 线性 常 微分 方程 边 值 
问题 的 算 例 。 精 确 的 计算 结果 ， 表 明 论 文 所 提出 算法 求解 这 类 问题 的 可 行 性 。 为 进一步 理论 研 
究 与 工程 应 用 提供 良好 条 件 。 
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Optimization Algorithm of Boundary Value Problem for Nonlinear 
Ordinary Differential Equations 


HOU Xiang-lin, QIAN Ying, WU Hai-tao 


(School of Science, Shenyang Jianzhu University, Shenyang 110168) 


Abstract: In this paper, an optimization method is applied to accurately solve a boundary value 
problem of ordinary differential equations. In view of the weakness of the shooting method, this 
study treats the initial conditions of the unknown parts as design variables and conducts complicated 
objective function based on the relationship between the given boundary values and design variables. 
And therefore, a new optimization algorithm to compute initial conditions of unknown parts is proposed, 
which is finally verified by employing Visual Basic 6.0 to design a universe program, providing typical 
examples and making comparison with the existent results. 
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